Interaction beween polarons and analogous effects in polarized Fermi gases 
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We consider an imbalanced mixture of two diflerent ultracold Fermi gases, which are strongly 
interacting. Calling spin-down the minority component and spin-up the majority component, the 
limit of small relative density x — ni/n^f is usually considered as a gas of non interacting polarons. 
This allows to calculate, in the expansion of the total energy of the system in powers of x, the terms 
proportional to x (corresponding to the binding energy of the polaron) and to x^^^ (corresponding to 
the kinetic energy of the polaron Fermi sea). We investigate in this paper terms physically due to an 
interaction between polarons and which are proportional to and x'^^^ . We find three such terms. 
A first one corresponds to the overlap between the clouds dressing two polarons. The two other 
ones are due to the modification of the single polaron binding energy caused by the non-zero density 
of polarons. The second term is due to the restriction of the polaron momentum by the Fermi sea 
formed by the other polarons. The last one results from the modification of the spin-up Fermi sea 
brought by the other polarons. The calculation of all these terms is made at the simplest level of a 
single particle-hole excitation. It is performed for all the possible interaction strengths within the 
stability range of the polaron. At unitarity the last two terms give a fairly weak contribution while 
the first one is strong and leads to a marked disagreement with Monte-Carlo results. The possible 
origins of this discrepancy are discussed. 

PACS numbers: 03.75.Ss, 05.30.Fk, , 67.85.Lm, 71.10.Ca 

I. INTRODUCTION 

Ultracold atomic gases have proved to be remarkably interesting systems due to the extreme simplicity of their 
effective interactions, which can indeed be described in many cases by the mere knowledge of the scattering lengths 
between the involved atoms. In particular fermionic gases made of two kinds of atoms, usually two hyperfine states 
of a same atomic element, have their interactions fully described by the scattering length a between different atoms. 
A large number of studies have been devoted to the case of balanced systems where the numbers of the two different 
atoms are the same [l]. More recently attention has also turned toward the case of imbalanced mixtures, where these 
atom numbers are unequal and which display very rich physics. They are analogous to systems found in other fields 
of physics, in particular superconductors in very high magnetic fields, but the ability to vary at will the ratio between 
the two atomic populations makes them especially convenient to explore experimentally. 

A particularly interesting situation is the limiting case where one atomic species is very dilute compared to the 
other one, so that its understanding reduces to the study of a single fermion in the Fermi sea of the other population. 
We follow the fairly standard convention of calling this isolated atom a "spin-down" atom while the Fermi sea is 
made of "spin-up" atoms and consider only the case where this interaction is attractive. This interaction between 
the spin-down and spin-up atoms "polarizes" the Fermi sea and the resulting complex object, a quasiparticle, is quite 
often called a "polaron". It is characterized by its binding energy £'{, = — /i^, where is the chemical potential of 
the spin-down atom, and also by its effective mass m* which gives the kinetic energy /2m* when this quasiparticle 
has momentum p. This polaron has been studied by various theoretical approaches [5]-[7], variational, diagrammatic 
and Monte-Carlo, which have been found in very good agreement with experiments [8l-fl0]. 

The study of the polaron has been the basic step in the theoretical investigation by Lobo et al [2] of the first order 
transition observed for trapped imbalanced Fermi gases by the Rice and the MIT groups [TT]. Their study has been 
performed at unitarity, but it has been extended later on to the whole BEC-BCS crossover [T^]. They have evaluated 
the total energy of a normal gas of polarons in a Fermi sea of atoms by taking into account not only the 
binding energy of these polarons, but also the kinetic energy of the polaron Fermi sea, proportional to x'^/^, where 
X — nj^/n^ is the spin-down concentration compared to the spin-up one. These two ingredients have proved essential 
to the analysis of experiments and the quite good agreement between experiments and theory. More generally the 
approach starting from the dilute limit for the spin-down population and leading for the total energy of the system 
to an expansion in powers for the spin-down density x has proved to be extremely fruitful. Our paper is in this line 
of reasoning. 

However Monte-Carlo calculations also performed by Lobo et al [2^ showed a small and progressive departure for 
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larger concentration x from the energy just given by the sum of the polaron chemical potential and the kinetic 
energy of the polaron Fermi sea. Such a departure is naturally expected since, when the polaron concentration is 
increased, the dilute regime, where each polaron can be considered as isolated, breaks down and one has to consider 
the contribution coming from the interaction between polarons. Once the properties of the isolated polaron is under 
control, this interaction between polarons is naturally the next problem of theoretical interest. But it is also of 
practical interest in order to have a better understanding of the agreement between theory and experiment. This is 
the basic problem we will consider in the present paper. We note that Nishida 17 as well as Patton and Sheehy [T3] 
have investigated the quite interesting question of a possible instability coming from polaron interaction and leading 
to the p-wave superfluidity of the polaron gas. However we will restrict ourselves in this paper to the calculation of 
the total interaction and we will not proceed to an analysis of the partial wave components of the effective interaction. 

More precisely this analysis in terms of polarons is valid only in the dilute regime of small x, strictly speaking in the 
limit X — > 0. In this regime we can meaningfully expand the total energy of the system in powers of x. The validity 
of this power expansion is only guaranteed for small x, although one may hope, as it happens quite often, that it 
extends actually to fairly larger values of x including for example the value Xc — 0.44 where Lobo et al have found 
a first order transition at unitarity. The binding energy of the polaron gives in the total energy a term proportional 
to X and the kinetic energy of the polaron Fermi sea is proportional to x^^^ . Hence we are interested in terms with 
next higher power. These will turn out to be terms proportional to x'^ and x'^^'^. It is worth noting that, although 
these exponents are indeed higher than the 5/3 arising for the polaron Fermi sea kinetic energy, the corresponding 
terms are in practice quite relevant for quantitative evaluations since, in order to neglect them, one would need to 
have typically x^^^ ^0.1 which implies extremely small values of x ^ 10^'^, that is much smaller than the typical 
values Xc for which the first order transition has been found at unitarity. 

The most obvious physical origin of the interaction between two polarons is the overlap of the spin-up clouds 
surrounding the spin-down fermions, that is the overlap between the polarized regions forming the polarons. This 
makes each of these clouds less "perfect" than for an isolated polaron. This reduces the binding energy of each polaron, 
which implies a repulsive interaction between these polarons. However, in contrast with the naive expectation, the 
result for the total energy is not proportional to x^, but rather to x'^/^ . This is directly due to the Fermi exclusion 
principle between polarons as we will see. Indeed, qualitatively, this makes it impossible for two polarons to be in the 
same state, and depresses their probability to be spatially nearby, which decreases their overlap. 

Apparently this argument leads to the surprising conclusion that there is no contribution proportional to x^. 
However this is clearly incorrect. Indeed such a contribution arises because, when we calculate the chemical potential 
of a single polaron and want to have our result for the total energy correct up to order x^, we have to take into 
account that the Fermi sea of the spin-up fermions is modified at first order in x by the presence of the spin-down 
fermions, which gives rise to a correction of order x to the chemical potential of the polarons in the dilute regime, 
and consequently a contribution of order x^ for the overall energy. Moreover, in order to obtain a proper expression 
of the polaron binding energy, we have also to take into account the modification of the spin-down properties, namely 
the existence of the polaron Fermi sea, which is enough at the order we are working. One could expect that this last 
correction gives another contribution of order x^ for the total energy, but it happens to be of order x^^'^, that is of 
the same order as the term resulting from the overlap of the polarization clouds. 

In the following we will calculate in turn these three above contributions. We will systematically stay at the simplest 
level of approximation, namely restricting the spin-up Fermi sea excitations to single particle-hole excitations. This 
is first for simplicity since the involved calculations are already not so simple. On the other hand this is also justified 
by the fact that this approximation has proved to be excellent for the calculation of the polaron binding energy and 
for its effective mass [SHS HI H] , for reasons which are quite well understood [7] . We may expect these same reasons 
to carry over into our problem, although we will not attempt to prove it. There is no difficulty in principle to go to 
the level of two particle-hole excitations or more within the framework we use, but in practice this would imply much 
more work. In the concluding part we will come back to the discussion of the power expansion and also compare our 
results to recent works by Mora and Chevy [T^ , and by Yu, Zollner and Pethick [T^] who have adressed basically the 
same problem as the one we consider here. 

II. OVERLAP OF POLARIZATION CLOUDS 

Our purpose in this section is to calculate the first correction to the dilute limit, corresponding to the term with 
lowest power in x related to the effect of the overlap between polarization clouds. We will proceed in the following 
way. We will first consider two polarons with momenta pi and p2 and calculate their total energy. Their interaction 
energy Z//(pi, P2) is the difference between this total energy and the sum of the energies of the isolated single polarons 
with momenta pi and P2 respectively. Then the contribution E^^^ to the total interaction energy coming from the 
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overlap is obtained by summing up over the free polaron Fermi sea, with radius pp: 

Elnl = l E "(Pl'P2) (1) 

P1;P2<PF 

with the factor 1/2 to avoid double-counting. Clearly this leaves out a number of effects which would lead to higher 
order corrections. For example the interaction energy of three polarons is not in general obtained by merely summing 
the three polaron-polaron interactions, there are specific three-body contributions. Similarly we are allowed to consider 
the free polaron Fermi sea only at lowest order, since in general the Fermi sea itself will be modified by interactions. 
However taking into account such an effect would again lead to higher order terms in the expansion of the energy in 
powers of the density x. 

In order to perform the above program we will extend our preceding work [7|, where we were dealing with a single 
spin-down interacting with a free spin-up Fermi sea, to the case where we have two spin-down fermions. In principle 
we would like to write and solve the full corresponding Schrodinger equation. Actually this full problem is much too 
difficult and we will restrict ourselves to a much simpler one by restricting as much as possible the Hilbert space. In 
addition to the two spin-down fermions we allow only the presence of a single particle-hole excitation in the spin-up 
Fermi sea. Restricting in such a way the Hilbert space is equivalent to perform a variational calculation. Naturally the 
first justification for such a reduction of the Hilbert space is simplicity. However we believe that we will nevertheless 
obtain a fairly good result for polaron interaction. Indeed it has been shown that, for the single polaron problem, 
allowing a single particle-hole excitation provides an excellent result [3j IH [7] for the polaron chemical potential. It is 
reasonable to believe that the reasons behind the success of this approximation [7] will also work for the case of two 
polarons. However we will not try to prove it because we would need to allow at least two particle-hole excitations 
which, although quite possible, is much more complicated than the calculation we will perform here. 

Accordingly we consider for the two polarons an eigenstate of the form; 

k^L^CqlO) (2) 

PlP2 piP2kq 

where |0) = rifc<fcF ""k k'^c) is the Fermi sea of f-spins, Ck and are annihilation and creation operators for f-spin 
atoms while bp and &^ are the corresponding operators for the 4--spin atoms. In the following we assume implicitly 
k > kp and q < kp. Since the operators b^p_^ and 6^^ anticommute, the coefficient Ofp^pj is naturally antisymmetric 
with respect to the exchange of pi and p2: 

'-^PlP2 ~ ~'^P2Pl (3) 

implying the Pauli exclusion principle. Similarly in the second term, by anticommuting the two operators and 
making the appropriate change of variable, we conclude that: 

Q^piP2kq = ^Ctp2+k— q,pi+q— k.kq (4) 

Hence we have for example in the first term of Eq.([2| two identical terms describing a spin-down fcrmion with 
wavevector pi together with a spin-down fcrmion with wavevector p2. To correct for this redundancy we should put 
a factor 1/2 in front of the eigenstate Eq.([2]). However since the normalisation of the wavefunction is unimportant in 
the following, it is simpler not to write this corrective coefficient. 

The Hamiltonian of our problem is the standard one channel Hamiltonian 1 , widely used to investigate the single 
polaron problem, namely: 

H = H, + V (5) 
= ^i?(p)&^6p + Eek4ck 

P k 

^ = 5 E ^kk'pp'4ck'6j,V 

kk'pp' 

where ek = k^/2m^ and E{p) = p^/2m^ are the kinetic energies (we take h = 1 throughout the paper). In the 
interaction energy term V, the Kronecker symbol insures momentum conservation. Explicitly (5kk'pp' = <5k+p-k'-p'- 
In contrast with the single polaron case, where it is convenient to set the volume equal to unity, it is here necessary 
to write explicitly the volume V in all the formulae. Indeed when we let this volume go to infinity, we will have 
physically two independent polarons with a total energy equal to the sum of the energies of each polaron. If they have 



4 



Pi and P2 as respective momenta, their energies are respectively fi +p\/2m* and /z +p\/2m* where fj, = m and m* 
are the chemical potential and the effective mass of a single polaron [SHZ], which depend naturally on the masses of 
the atoms, on the density n = kp/Qn^ of the spin-up atoms and on the scattering length a. On the other hand if the 
volume V is large but finite we expect, in addition to the preceding contributions to the total energy, an interaction 
energy which will scale as l/{kpV), since the typical volume of the polarization cloud of spin-up atoms surrounding 
a spin-down one is of order 1/kp, which has to be compared to the sample volume V. Hence we have to consider a 
finite volume if we want, for the interaction energy we are interested in, to recover a nonzero result. 

This can be done in two equivalent procedures. A standard one is to introduce a factor 1/V in the right-hand side of 
the interaction term in Eq.([5]). In this case the dimension of g changes and, instead of being an energy, it becomes an 
energy multiplied by a volume. A completely equivalent, and for our purpose slightly more convenient, way is to keep 
Eq. ([5]) unchanged so that g has still the dimension of an energy. In this case the relation between the scattering length 
a and the interaction g has to be modified. For a unit volume this relation reads mr/(27ra) — g~^ + X)'^'^ 2mr/k^, 
where rur = m^mi / {m^ -\- m^) is the reduced mass and k^ is an upper cut-off in the interaction energy necessary to 
avoid ultraviolet divergences. As usual we will let this cut-off go to infinity while the coupling constant g goes to zero, 
keeping the scattering length finite in the preceding relation. Physically this corresponds to take the limit of a very 
short ranged interaction potential. However in the general case the dimensions of the preceding formula are incorrect 
and we have to write instead: 

2™ ~ g ^ k^ ^ ' 

where one sees that all the terms have the same dimension. This equation is clearly equivalent to the equation 
obtained by the first procedure (which is merely obtained by dividing by V and making the appropriate change of 
notation). It may also be checked in the situation where = 00, in which case we have merely a scattering potential 
located at the center of a box of volume V. This equation is convenient because nothing is changed compared to the 
case where the volume is unity, except that we have to make the single simple replacement a — !■ a/V. 



A. Projected Schrodinger equation 

Now, we write the Schrodinger equation H\iIj) = E\ip) and project it on the subspace corresponding to Eq.([2|, 
i.e. onto the full Fermi sea and the Fermi sea with in addition a single particle-hole pair. This yields a set of two 
coupled equations. We will not write them for the general case, but rather take advantage of simplifications which 
were already arising in the case of a single polaron j^- First some terms, like the Hartree term, disappear in the limit 
of vanishing interaction strength g — >■ which we have to take, as we have explained above. Naturally we will not 
write these terms. Second the terms resulting from the scattering of an spin-up particle (meaning with wavevector 
larger than kp) display ultraviolet divergences analogous to the one present in the last term in the right-hand side of 
Eq.Q, linked to the fact that the wavevector can go to infinity. These will be precisely cured by making use of Eq.(|6|. 
On the other hand terms physically related to the scattering of holes (which have wavevectors smaller than kp) do 
naturally not present these ultraviolet divergences. Hence in the above limit of kc — ^ 00,17 ~^ they are negligible 
compared to the preceding ones. Hence we will take only into account terms corresponding to the scattering of spin-up 
particles and will not write those corresponding to the scattering of holes, because in the end these last ones turn out 
indeed to be vanishing as can be checked by keeping them all the way long. With these simplifications we end up 
with: 

- 5"^apiP2 + E{-p2) - E)]^^ ("pipakq - ap2Pikq) (7) 

kq 

and 

- 5"^apiP2kq -E'pipakq = '^PlP2 - ap2+k-q,pi+q-k + ^ ("pipsKq ^ "ps+k-q.pi+q-k.Kq) (8) 

K 

where we have set: 

4?P2kq = ek - eq + i?(pi + q - k) + E{V2) - E (9) 

(we expect naturally E <Q). We note that the antisymmetry properties Eq. ([s]) and Eq. Q are automatically satisfied 
by these equations. This has been obtained by making use, when necessary, of the anticommutation property of the 
two operators. 
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Then we manipulate and simplify these equations in much the same way as we have done in the case of a single 
polaron [7]. Indeed the complication present here is the dependence of the wavefunction on the spin-down variables 
Pi and P2. However the key part of the manipulation comes from the dependence on k, in particular for large k, and 
this aspect is essentially unchanged. Indeed it is clear from Eq.([7]) that ap^pjkq ~ since the summation over 

the hole variable q plays an unimportant role. Dividing Eq.(j8| by E^p^^^y^^ and summing over k leads to: 

k \ K / \ k ^pipakq / 

~ 51/ T^(T) I '^P2+k-q,Pi+q-k + ap2+k-q,pi+q-k,Kq 1 

k ^pip2kq \ K / 

where we have introduced Eq.([6]). The dominant contribution g^^ Opipakq cancels out. Then in the last parenthesis 
the first term ap2+k-q,pi+q-k is negligible compared to the last one Op^+t-q.pi+q-k.Kq ^ 9^^ in the limit g — > 0. 
We are left with: 

(\ ^ 1 Vrrir ^ 2mr\ sr-^ 

^ + ^ - L ^ L "P^P^Kq (11) 

k -'^pipakq / K 

= .9 CtpiP2 + 5^ 3(1) 5^ CKp2+k-q,pi+q-k,Kq 

k ^pip2kq K 

where in the factor of Q;pjp2, in the first term of the right-hand side, we have kept only the dominant term —5^^, the 
other ones being negligible in the limit g — > 0. Setting: 

^q(Pi,P2) = sXl^PiP^Kq (12) 

K 

which is finite in the limit g — > 0, and with the notation: 

[T,(pi,p,)] ^^^_ + _-^_ (13) 

k ^pip2kq 

we end up with: 

_ Fci{Pi,P2) _ ^q(P2 + k q,pi+q k) 

''^^p^"rq(Pi,P2) E^'^ k 

^ k>kF -^pip2kq 

If we look at the volume dependence of the various terms, we see from Eq.( |13[ ), after making the conversion from 
summation over k to integration over k by — > V/{2Tr)^ J dk, that [rq(pi, p2)]^^ is proportional to volume V. 

It is convenient to remark that this equation relates terms in the wave function corresponding to a fixed value of 
P = Pi -|-p2, which is natural since the total momentum is conserved in the scattering of the two spin-down particles. 
Setting pi = P/2 + s and p2 = P/2 — s, with ap(s) = apip2 and Fqp(s) = i^q(pi,p2), and similar notations for T 
and E, we can rewrite Eq.([T4|: 

- (.\ ■^qP(^) fqp(k-s-q) 

' k>kF -f^Pskq 

For fixed P and q this is just a matrix relation between the vectors ap(s) and Fqp(s). Making the change 
k = s + t + q, we can write: 

ap(s) = V^Aqp(s,t)Fqp(t) (16) 
t 

where the matrix elements of A are given by: 

A '^s,t 1 g(|s + t + q|-fcf) 

^ ' ^Ps,s+t+q,q 
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and we have set: 



7qp(s) 



V 



tqp(s) 



Since: 



41+t + q,q = ^s+t+q " + ^ s) + E - t) - E = 4M+t + q.q 



the A matrix is symmetric: 

With our new notations Eq.([7]) reads: 

Eap(s) = 



v4qp(s,t) = Aqp(t,s) 



'-E{F)+2E{s) 



(18) 



(19) 



(20) 



(21) 



where the two first terms (l/2)i?(P) + 2E{s) — E{pi) + E{p2) represent the kinetic energy of the two spin-down 
particles written as the sum of the kinetic energy associated to the center of mass motion and the kinetic energy of 
the relative motion. Inverting Eq.(16) and carrying it into Eq.(21| we obtain: 



E - -EiP) - 2Eis) 



qt 

= ^ E (^qP *) + V (-^' -*)) (*) 



(22) 



where we have used Eq.([3]) in the last step. We will now make use of the explicit expression Eq.(17) for Aqp(s,t). 



Infinite volume 



We consider first the infinite volume limit where the second term in the right-hand side of Eq.(17) is zero. In this 
case A is diagonal and its inverse is: 



^qp(S't) = iqP(s) <5s,t 



Carrying this expression into Eq.(22) and going to integration over q we have: 

1 



E - ^EiP) - 2E{s) 



ap{s) 



(27r)= 



dq (tqp(s) + tqp(-s)) ap(s) 



or with our original notations (with iq(pi,P2) — iqp(s)): 

[E - E{pi) - E{p2)] ap^p^ = y"rfq (iq(Pi,P2) +^q(P2,Pi)) apip2 



(23) 



(24) 



(25) 



The volume has disappeared from the equation and at first this result looks like what is expected for the energy 
of two non-interacting polarons. Indeed let us consider the first term in the right-hand side. It gives the chemical 
potential of the first polaron pi. Indeed for this polaron the second polaron p2 is irrelevant, and the relevant energy 
for this first polaron is E~E{p2), that is the total energy E from which the kinetic energy E{p2) of the second polaron 
has been removed. Similarly one can check that this is also this difference E — E{p2) which is entering tq(pi,P2) 
through Eq.Q and Eq.(13). With the change in notations E — E{p2) E, one recovers for this polaron: 



E = E{pi) 



dk 



(27r)^mr 



Ek - eq + E{pi +q-k)-E 



2mr 

'>' 



(26) 



which is exactly [21 HJ [7j the implicit equation for E which gives the single polaron chemical potential together with 
its effective mass. Finally we can see in the same way the second term in the right-hand side of Eq.(25l as giving 
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the chemical potential of the isolated second polaron p2 since, by exchanging the variables pi and P2 and taking the 
antisymmetry Eq.([3]) of Ofp^pj into account, it becomes the first term of the equation and we can repeat the above 
analysis. 

Actually the above argumentation is clearly not completely satisfactory since the analysis of the first polaron energy 
implies that E = /i + p^/2to| +p2/2mi while the second one corresponds rather to E — n +p\/2m*^ +p\/2m^. Since 
in general |pi| 7^ IP2I, these two values are inconsistent. We might perhaps argue that, since we are anyway interested 
in small values of |pi| and |p2|, we should neglect kinetic energies and only the zcroth order result E — jiis relevant. 
However even in this case we see that the sum of the two contributions in the right-hand side of Eq.(25) gives 2/i 
instead of /i, so that this equation is actually not satisfied. 

The reasons for these problems are physically quite clear. Since we have taken the very strong restriction Eq.([2| for 
the Hilbert space, we do not even allow the situation where each one of the polarons has its own single particle-hole 
dressing. Although it is not in general meaningful to attach the particle-hole pair to a specific spin-down particle, we 
may roughly say that when one spin-down particle is dressed by the particle-hole pair, giving it the energy jj, (for zero 
kinetic energy), the other spin-down particle is naked and its energy is just its kinetic energy. In this perspective the 
expressions written in the above paragraph for the total energy E become quite meaningful. 

A tempting way to solve the preceding inconsistencies is to enlarge our Hilbert space by allowing the existence of 
two particle-hole pairs. In this case we would have the possibility to have each spin-down particle dressed by one 
particle-hole pair, allowing a much more satisfactory physical description of the two polarons. However this would 
not solve all the consistency problems since we would have the possibility that one down-particle is dressed by two 
particle-hole pairs while the other down-particle is naked. This would be a situation analogous to the one we are faced 
with, although it would clearly be much more complicated since we should treat the single polaron problem at the 
level of two particle-hole pairs. It is nevertheless possible to check that, if we allow each spin-down particle to have a 
single particle-hole pair dressing, the above inconsistencies disappear. Indeed we can consider the case where the two 
spin-down particles belong to different systems, each spin-down particle polarizing a different Fermi sea. Hence one 
has two non-interacting systems. Physically this is not so different from the case of two polarons in an infinite volume, 
since the overlap between the two polarons will also be zero in this last case. The state describing this situation is 
merely the tensorial product of the two polarons (each one having its single particle- hole pair dressing) . Hence there 
is the possibility to have two particle-hole pairs. One can then write the corresponding Schrodinger equation, which 
decouples naturally into two Schrodinger equations, describing the two independent systems. Hence one can follow 
how the existence of the additional term with two particle-hole pairs solves the above inconsistencies and allows the 
decoupling of the equations, leading to the obvious result that the total energy is the sum of the energies of each 
polaron, each one being obtained by solving Eq. ( 26 1 . Since writing all this explicitly does not make any problem, but 
is fairly lengthy and burdensome, the reading being uninspiring, with a result which is completely obvious physically, 
we will not display the corresponding equations and stay at the level of this physical discussion. Our final conclusion 
is that, since at this stage everything reduces to the calculation of the energy of a single polaron, the energy E which 
comes in the expression of tq(pi,P2), through Eq.Q and Eq.(13l, is the single polaron energy, i.e. the polaron 
chemical potential /i when we work at zeroth order in kinetic energies. 

As a final remark we note that all the problems discussed above are hidden in the case of the weak coupling limit, 
where iq(pi,P2) = Ina/nir and each term in the right-hand side of Eq.(25) gives merely a (mean-field) contribution 



C. Polaron interaction 



Let us now consider the case of a finite volume V. In this case the second term in Eq.(|17[) is no longer zero. However, 
because of the factor 1/V, all the corresponding matrix elements are very small in the limit of a very large volume 
compared to the matrix elements arising from the first term. This makes the inversion of the A matrix very easy, all 
the more since the first term corresponds to a diagonal matrix. One obtains: 



Aq^(s,t) = fqp(s)4.t + y tqp(s) ^^'^ 



E, 



^ fqp (t) = tqp (S) 5,,t + ^ 5qP (s, t) 



(1) 

Ps,s+t+q,q 



where is the Heaviside function. Upon substitution in Eq.(|22|) we find 
1 



E - -E{V) - 2E{s) 



1 



aplsl 



(27r)3 
1 



dq [tqp(s) + <qp(-s)] ap(s) 



+ V 5 (2^) 



lyi J dcL [Bqp(s,t) +Bqp(-s,-t)]ap(t) 



(27) 



(28) 
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Making use again of Q;p(— t) = — ap(t) the last term in the right-hand side may also be rewritten as: 

/ [Sqp(s,t)-i?qp(s,-t)-i?qp(-S,t)+i3qp(-S,-t)]ap(t) 



(29) 



Let us now consider what would happen if we had two spin-down particles interacting directly through a general 
non-local interaction V and being ruled accordingly by the Hamiltonian: 



p KK'Q 

where, because hermiticity and antisymmetry under fermion exchange, we must have: 

Vq(k,k') = y4(K',K) = -^q(-k,k') = -yQ(K,-K') 



When one writes the corresponding Schrodinger equation, one finds instead of the last term of Eq.(28|: 

^Y.^VAs,t)-VF{s,-t)]aF{t) = ^ ^yp(s,t)ap(t) 



(30) 



(31) 



(32) 



where we have made use of the antisymmetry relations Eq.(31 1. 



When we compare Eq.(32) with Eq.(29), we see that we can interpret Eq.(28) as the Schrodinger equation for the 



two spin-down particles interacting through an effective potential given by: 

^P(s,t) = ^^^ / '^'^ [5qP(s,t)-Bqp(s,-t)-Bqp(-S,t)+Bqp(-S,-t) 



(33) 



We note that this expression satisfies as it should the requirements Eq.(31 ). Naturally this effective interaction is due 
physically to the fact that the two polarons have to share some part of their polarization clouds. 

It is worth noting that our approach is quite general. Restricting ourselves to a single particle-hole pair leads to the 
explicit expression Eq.(27) for B. However, if we extend our Hilbert space and consider any number of particle-hole 
pairs, Eq.([7]) is still valid. And by eliminating the wavefunction components corresponding to more than a single 
particle- hole pair, we will end up with an equation similar to Eq.(16). The explicit expression of matrix A will be 
more complex, but otherwise we can proceed in exactly the same way as we have done above to identify the effective 
interaction, and consequently obtain the interaction energy as we will do just below. Hence extension to taking for 
example into account two particle-hole pairs can proceed in this way, with the same framework. 

If we calculate now the mean value of the interaction energy in Eq. ( 30 ) for a state where the two spin-down particles 
occupy plane waves pi and p2, namely the state |tiac), we obtain easily: 



V 



¥p(s,s) 



(34) 



with again the notations P = pi -I- p2 and s = (pi — P2)/2. Making use of Eq.(33) this leads for the two polarons to 
an interaction energy: 

1 1 f 

Z^(Pl,P2) y dq [Bqp(s,s) - Bqp(s, -S) - Bqp(-S,s) -|-Bqp(-S, -S)] (35) 

Making use of Eq.([T]) we obtain finally: 

^int=^ E "(Pl'P2) = ^ E ^ /^^q [SqP(s,s)-Bqp(s,-s)-i?qp(-S,s)+Sqp(-S,-s)] (36) 



^1,P2<PF 



Pi,Pi<Pf 



We could as well calculate directly the average of the interaction energy Hint given by Eq. ( 30 1 over the spin-down 
Fermi sea |0) = np<pj7 ^p l^'^c). We would obtain a direct term and an exchange term, leading to: 



(1) 



^ ^ [Vi,(s,s)-l/p(s,-s)] = i 



2V 



Pi,P2<Pf 



(37) 



Pi,P2<Pf 



which is identical to Eq. ( 36 ) from Eq. ( 33 ) . 
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D. Total interaction energy 



Let us now proceed to evaluate explicitly this interaction energy. We note first that, by going from summation 
to integration through ^ — V/{2Tr)^ J dpi, we find that E^^^ is as expected an extensive quantity, proportional 
to volume. Having settled this point, we will take in the following V = 1 for simplicity. Then we notice that the 
dependence of E^^^ on will be faster than the standard dependence in n| ~ expected for a standard interaction. 



This is due to the fact that Vb(0, 0) — 0, so the naive answer E' 



(1) 



?Vb(0,0) does not hold. This cancellation, seen 



explicitly in Eq.(33), results automatically from the antisymmetry properties Eq.(31| and so it comes directly from 



the fact that polarons are fermions, just as the spin-down particles themselves. 

We will naturally use the fact that is small compared to kp, and for consistency we have to make the calculation 
to lowest order in pp. However, in view of the above mentioned cancellation, we have to do it carefully. Looking 
at the explicit expression Eq.(27) of Bqp(s,t), we see that if we set s = t = in the Heaviside function we obtain 
9{q ~ kp), so that q must be larger than kp. However by definition q < kp. Hence the result will be zero due to this 
Heaviside function. This means that, although |s|,|t| < pp we have to keep carefully nonzero s and t in this Heaviside 
function. On the other hand we may perfectly set s = t = in the other factors since they are not sensitive to this 
approximation and the result is nonzero. Similarly |P| < 2pp allows us to set P = since the resulting values of the 
different factors are nonzero. A further simplification comes from the fact that i?qp(s,— s) — i?qp(— s, s) = since 
both are proportional to 9{q — kp) which leads to a zero result as explained above. Finally exchanging pi and P2 
changes s into — s, hence the last term in Eq.(36l gives the same result as the first one. We are left with: 



E.. 



(1) 



[ikpo(0)]2 



1 



i2n) 



dp I 



Pi<Pf 



P2 



Pi<Pf 



dq6'(|pi -P2 +q| - kp) 



(38) 



q<kf 



where we have used i?oo^q.q — \fJ-\ from Eq.(|9|), since we have explained at the end of section 



HB 



that we have to take 
p from the Heaviside 



E = II. We have also replaced iqo(O) by tkFo(O), because to zeroth order in pp we have q 
function. 

Hence we are left to evaluate the high dimensional volume corresponding to the triple integration. The q integration 
is easily performed. Indeed q must be inside a sphere of radius kp centered at the origin, but from the Heaviside 
function outside a sphere of radius kp centered at 2s = pi — p2. Thus the result is the volume inkp/S of the whole 
sphere minus the volume of their intersection. Since the volume of the intersection of two spheres of radius R with 
their center distance being r is given by Vn = (7r/3)(4i?^ — 3R^r + r^/A), this gives: 



f 27r 

/ dq 0{\pi - P2 + q| - kp) = — (3s/c| - s^) ~ 2nsk% 

Jq<kF 



(39) 



where we have used s ^ fcj;^ in the last step. Setting pi — pi?x and p2 = —ppy, we are left, apart from a factor p'p, 
with the calculation of: 



dx 



X<1 



dy |x- 



y<l 



(40) 



for which it is more convenient to take the new variables u = x + y and v = (x — y)/2, with dxdy = dudv. Since it 
is easy to see that, at fixed u, vector v (with its origin at the middle of vector u) has to run inside the intersection of 
the two spheres of radius 1, with distance u between centers, the v integration gives again Vd with R — 1 and r = u, 
i.e. (7r/3)(4- 3w + uV4)- Hence: 



[ dx [ dy |x + y| = ^47r /" 
Jx<l Jy<l "J Jo 



duu^{A - 3u H ) 



647r^ 
~35" 



(41) 



This leads us finally to: 



^int - 



1 647r^ 



^1 (27r)9 35 



35 V47r4 



1/3 



Kp 



(42) 



where tkp = ikFo(O). Hence we obtain an interaction energy which scales as tJJ^, and not as rt^ as discussed above. 
One can express tkFo(O) analytically in terms of kp, a, m^, and but there is no point to write here the 



T 
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FIG. 1: (Color online) e';^^^ / in^Epx' as a function of 1/kpa in the case of equal masses = 7114,. The dashed line is the 
corresponding expansion in the weak coupling limit. 

somewhat comphcated formula. In the weak coupling limit a — > 0_, it reduces to 2'Ka/mr. Since in this case we have 
also the mean field value = 2'Kn^\a\/mr^ we obtain explicitly: 



(1) 



kl 



9 , 

pp\a\ 

rrir 



35tt 



(43) 



Let us finally note that this interaction energy 
polarons. This is already seen at the level of Eq.(38|) 



is positive which corresponds to a net repulsion between two 
This result is physically reasonable: two polarons have to 



compete to make up their spin-up clouds, so that their dressing is not as optimal as the cloud of a fully isolated 
polaron. 

In order to compare the above interaction energy to other contributions to the system total energy, it is convenient 
to write the dimensionless ratio e[^^ /n^Ep, where Ep = /i^ = k'^p/2m^ and to introduce the ratio x = n^/n-^. This 
leads to: 



E 



(1) 



3 [m^kptk^ 



n^Ep 35 TT^ 



(44) 



where we have set = pEp. In the weak coupling limit, this result becomes 18m.-fkp\a\ x'^^'^ /{35nmr). For the case 
of equal masses = — 2mr, the coefficient in front of a;^/-^ in the result Eq.(44) is plotted in Fig|lj We have 
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limited the plot range to l/kpa < 1, since we know [51 [TB] that, for l/kpa = 0.88, the polaron becomes unstable with 
respect to the formation of a (dressed) bound state between the spin-down and a spin-up fermion, so that beyond 
this point the whole physical picture breaks down. 



III. CORRECTIONS TO THE SINGLE POLARON CHEMICAL POTENTIAL 



Let us now come to the corrections we have to make for the single polaron chemical potential, due to the presence of 
a non zero density of spin-down atoms. As we have already indicated in the introduction, they have two origins. The 
first one is the modification to first order of the properties of the spin-up Fermi sea, which implies a corresponding 
change of the spin-up propagator coming in the calculation. The second one is merely the existence of the polaron, 
or equivalently spin-down Fermi sea. We begin by considering this last effect since it is by far the simpler one. 



A. Correction due to the polaron Fermi sea 



When one calculates the chemical potential of a polaron, one assumes that, when a particle-hole pair is created, 
the spin-down atom may have any recoil to compensate for the momentum of this particle-hole pair. However in the 
presence of a Fermi sea of spin-down particles this is no longer allowed since Pauli exclusion forbids all the momenta 
inside this Fermi sea. Hence in Eq.(26) which gives the chemical potential we are looking for (provided we set the 
momentum pi = 0), we have to put in the sum over k a further restriction |k — q| > pp. The resulting small 
modification of the integral is: 



Sliq) 



I 



dk 



0{\k 



Pf) 



ek-eq + E{q-k) + \fi\ 
9{pp ~\k~q\) 



dk 



Ck - fq + E{q - k) 



(45) 



dk 



ck - fq + £^(q - k) 



where we have naturally set E — —\fj,\ since we are looking for the lowest order correction. We have by definition 
|k| > kp and |q| < kp- Hence 0{pp — |k — q|) implies that |k| ~ |q| ~ kp since pp is small. To lowest order the 
denominator is just equal to |^| and we are left with the calculation of / dk9{pp — |k — q|) which is just the volume 
of the part of a small sphere of radius pp^ centered at q and which is outside the sphere of radius kp centered at the 
origin. Since pp is small we may consider the relevant part of this last sphere to be a plan, and we have just to find 
the volume of the pp sphere which is beyond a plan with nearest distance to the center of this sphere equal to kp ~ q, 
with q = |q|. This volume is easily found to be {Tr/3)[2pp + kp — q][pp — {kp — q)]^, with naturally the condition 
< kp ~ q < Pp, leading to: 



SI{q) = ~—-i2pp + kp 



q)iPF + q- kpf 9{pp + q~kp) 



(46) 



This modification (5/(q) induces a shift 5\iJ,\ in the chemical potential. Differentiating Eq.(26) with respect to E we 
find: 



dq«(q) 



Z (27r)6 

where we have used the fact that 5/(q) is nonzero only for q ~ kp. We have set: 



(47) 



(27r) 



dqfqo(0)]2 / dk 



E, 



(1) 

OOkq 



in which the k integration can be performed analytically, but the q integration has to be done numerically, 
integral in the right-hand side of Eq.(47) is easily evaluated. This leads to: 



1 
Z 



(48) 



The 



(49) 



The result is negative as expected, decreasing the binding of the spin-down in the spin-up Fermi sea, since the phase 
space for creation of a particle-hole pair has been decreased. 
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(2) 

In order to obtain the corresponding contribution E^^l to the interaction, we have to make use of the definition 
of the chemical potential fi^ = d£/dn^, where £ is the total energy of the system. Since we find that the above 
spin-down shift in chemical potential is proportional to p% ~ '^•4^^'^) corresponding term in the interaction energy 
is: 



(50) 



where is given by Eq.( |49p . This contribution to the interaction energy is again a positive. 

We note that this term is proportional to pp. This is exactly the same dependence as the one we have found for the 
interaction between polarons (and in the same way it increases the total energy). Hence this contribution should not 



be overlooked. One could have naively expected a p|. dependence for 5\pt\ in Eq.(49), corresponding to the volume 
of the polaron Fermi sea. The additional pp factor comes from the fact that only a thin shell, of thickness pp, is 
involved in the hole phase space, corresponding to the variable q, and not the complete hole phase space. 

In the weak coupling limit a — >■ 0_ we have seen that ikp = 2iTa/mr and similarly iqo(O) = 2'iTa/mr, so the bracket 
in the right-hand side of Eq.(49) reduces to unity and we obtain for this contribution to the total energy: 



(2) 



nr 



-— PF|a|^ — 
IIJtt rrir 



(51) 



We see that it has exactly the same form as the interaction energy Eq.(43), with a coefficient which is not so different. 
Hence this effect ca n no t be omitted. 

(2) 

Finally, from Eq.(49|, the corresponding ratio E>^1 /n-^Ep is given by: 



E 



(2) 



3 1 [m^kpt]f^^Y 



.7/3 



n^Ep 112 tt" Z p 
The coefficient in front of x^/^ in Eq. ( 52 ) is plotted in Figj2] for equal masses 



(52) 



B. Correction due to the modification of the spin-up Fermi sea 



We consider now how the spin-up Fermi sea is modified by the spin-down population to first order in and how 
consequently the spin-down chemical potential p.^ is changed to first order in n^. This leads to a contribution to 
the total energy proportional to n^, that is exactly the dependence expected from an interaction between spin-down 
particles. The first stage is to find the spin-up self-energy. This is just a Hartree-like term coming from the spin-down 
population. The only difference is that, since the g goes to zero, we should not write merely a single interaction, but 
rather sum over repeated interactions, which leads to sum up ladder diagrams (see Figjs]). This gives: 



5]^(k, cj) = > 



k' 



duj' 
2i7r 



G;(k',c^')ro(k + k',c^ + c^') 



(53) 



In agreement with our above notations we have used again the convention = J (ik/(27r)'^. The lo' integration 
runs over the imaginary frequency axis. Here ro(k,a;) is the vertex used already in Ref.[4J to calculate actually the 
polaron chemical potential (it contains the bare propagator Gq^). Since the factor in Eq.([53|) will give rise to a 
factor n^, we have naturally to evaluate the sum of the ladder diagrams to zeroth order in n^, which is indicated by 
the subscript 0. Exphcitly: 



k 



2mr 
~p ^ o _L ,,. _L ,,0 



ffc - fit) 



n + fi-f + p^l-ek- E(k + K) 



(54) 



where the superscript indicates that the chemical potential has to be evaluated to zeroth order in n^. 
Tq{K,Q) is closely related to T'q(pi,p2) defined above in Eq.(13). 



Naturally 



Since we do not consider the possibility of a bound state between an up and a down particle, the only singularities 
of Tq{K, 17) occur only on the positive real axis (one can check that, for 17 < 0, the denominator in the integral in 
Eq. ( [54| is always negative) which correspond physically to the continuous spectrum of the scattering states of the up 
and a down particle. We calculate I]^(k, w) from Eq.(|53|) for imaginary frequency ui and then continue analytically 
the result for any frequency. In order to avoid the singularities of Fq we close the lu' path by a semi-circle at infinity 
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FIG. 2: (Color online) E^^^ /{n^Epx'^^'^) as a function of 1/kpa in the case of equal masses = mi. The dashed line is the 
corresponding expansion in the weak coupling limit. 



in the Re a;' < half-plane, so that the integration contour encloses this half-plane. The only contribution comes 
from the pole of G^(k',a;'). The down particle number is given by: 



E / S;G,(k',.') 



2i7r 



(55) 



where the integration contour is the same as above, enclosing again the Rew' < half-plane. In the case of a 
single spin-down particle, where is vanishingly small, the ground state corresponds to a pole of G^(k',a;') at zero 
momentum k' = and zero frequency w' = 0. In the case of a small number of spin-down particles, the poles will be 



similarly in the vicinity of k' = and w' = 0. Hence we can set k' = and w' = in ro(k -I- k', w -I- cj') in Eq.(53 ), so 
that we are left with an integral which is just Eq.(55 1. This leads us to: 



St(k,w) = n;ro(k,a;) 



(56) 



In particular in the weak coupling limit a — >■ 0_, where from Eq.(54| ro(i^; = 2T:a/mr, we find 



(57) 



which is just the expected mean-field result. 
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FIG. 3: Ladder expansion for the spin-up self-energy Eq.(53l. 



Having found the self-energy, we have the full spin-up propagator 

1 



Gt(k,c.) 



Got(k,a;) -I- i;^(k,a;)G5^(k,w) = Got(k, w) -f (5G'^(k, w) 



(58) 



where Go^(k,aj) — [uj ~ + fJ,^]~^ is the free spin-up propagator. Here we have in the last step expanded the result 
to first order in n^, consistent with our low spin-down density approximation. With this result we may, in our grand- 
canonical framework (which is just the T — limit of the standard finite T formalism), calculate the change Sn^ of 
spin-up particle density resulting from the spin-down particle density n^, at fixed /i^. It is given by: 



E 

k 



did 
2iTT 



(5Gt(k,w) 



E 

k 



I]t(k,u;)GL(k, 



dx 



E 



du) I]-|-(k,(u;) 
2i7r cj — Efc + a 



(59) 



Since, from Eq.(56), S^(k, cj) is analytical for Rew < 0, it is convenient to close again the contour around this negative 
frequency half-plane. The only contribution comes from the pole at a; — a;, which gives a factor 6{x — ek) since 
we want this pole to be at negative frequency. This leads to: 



6n^ = J2 



d(jj 



J(efe-/i^)E^(k,0) 



(60) 



where 9I]-|-(k, e/j — pL^)/duj is a short-hand for the partial derivative 9S^(k, w)/9cj taken for lo = — fJ-^- 

Interestingly, as pointed out in [T7_, Sn-^ has a simple physical interpretation for — >■ 0. It is merely the number of 
spin-up particles in the cloud forming a single polaron. This quantity, or equivalently v = dn^ / dn^ at fixed ii-^ which 
is obtained from Eq.(60) and Eq.([56|), satisfies a simple thermodynamic identity |T7]. We have ly = — where 
the derivative is taken at fixed n^. Taking — corresponds to the case of a single polaron and in this case in the 
simple framework used in Ref. Hi is given by: 



= E ^(^f ^ ej<')ro(K, €K - Mr) 



(61) 



K 



Taking the explicit derivative of this formula with respect to ^-f gives a result in full agreement with the thermodynamic 
identity. However, since from Eq.(54), ro(K, fl) has also a dependence on /i^ (as well as on /i^), we should take also 
these dependences into account in the calculation of the derivative. As a result the thermodynamic identity is not 



satisfied. This is not so surprising since Eq.(60) and Eq.(61) are only an approximate results. Nevertheless we have 
found numerically that the disagreement with the thermodynamic identity stays almost always quite small, whatever 
the value of the parameter l/kpa- The difference takes sizeable values only when one approaches the transition point 
l/kpa ~ 0.88 where the polaron becomes unstable with respect to the formation of a molecular bound state [51 118j. 

We consider now how the calculation of /i^ is modified when we take into account that Gq^ has to be replaced by 
Gq^ + (5G^. Taking again the framework used in Ref.[3], we have: 



Mi 



S,(O,0) = 5] / |^Gt(K,f7)r(K,17) 



(62) 
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Compared to the calculation in the single polaron limit — >■ which leads to Eq.(61), we have now to take into 
account that, in this formula, is no longer Gq^, but = Go-f + 5G^ given by Eq.(58). Similarly T is no longer 
Fq, but r = Fq + (5F, since in summing up the ladder diagrams, we have now to take into account that the spin-up 
propagator G^ is no longer Go-f. Since in this summation the contribution coming from the propagators between two 
interactions is basically a convolution of the spin-up and spin-down propagators, the modification will come from the 
quantity: 



E 

k 



2in 



GoiiK--k,n-uj)dG^ik,u;) 



E 

k 



duj 

2iTT 



Go^(K -k,n- u;)E^(k,w)G2^(k,u;) 



(63) 



where Go^(k, w) = [w + /ijj — E{k.)] ^ is, consistently with our above framework, the down propagator with the down 
self-energy set to zero, but with the proper polaron chemical potential. Since I]^(k,w) is analytical for Rew < 0, it 
is again convenient to close the contour around the negative frequency half-plane. Contributions will only come from 
the poles of Gq^ and Gq^ located in this domain. Using for the double pole due to Gq^ the same convenient trick as 



in Eq.(59), we obtain: 



E^(k,r! + ^J -S(K-k)) 



{n + ^i^ + fil- €k - E{K - k))2 dx^n + nl-€k- E{K - k) 



E 



S^(k, £k - x)9{x - ek) 



(64) 



One can check on this expression that (5(F ^)(K, ft), and accordingly F(K, il), is as expected analytical for Refi < 0. 
Indeed the apparent (double) pole at f2 = efe+£'(K— k) — /i^ — /i" does not exist since one can check that the numerators 
combine exactly to leave an analytical function at this point. This property allows to close again in Eq.(62 ) the contour 



around the negative frequency half-plane. From Eq.( |58[ ) this leaves us to evaluate the contribution coming from the 
poles of Gq^ and SG^, this last one being again handled as in Eq.(59). 



In obtaining the first order variation S^i of /i^ from Eq.(62 1, one last point has to be taken into account. Since there 



is a variation of /i^ and that, in the zeroth order contribution Eq.(61 1, Fq itself depends on fi^ as seen in Eq.(54|, we 
have to take into account this variation. In other words we have to keep in mind that Eq.(62) is an implicit equation 
for /i^. This introduces again the coefficient Z which has been introduced in Eq.(48l in the preceding subsection: 



K 



kK 



\E 



(1)1 

kKJ 



(65) 



where we have used the shortened version i?^K = -^ookK '^^ '^^^ notation Eq.([9|. This quantity is actually just 
1 — d'S^/dtd, that is the inverse of the residu e of G^ at its pole for zero momentum and energy. 

Making use of all these ingredients in Eq.(62), together with 5T = — Fq(5(F^1), we obtain for the variation Sfi^ of 
fj.^ to the first order in n^: 



Z 



^ e{fi-f - eA')Fg(K, ex - MT)ro(k, ex 

kK 



E{K - k) - Ait + m") d 



[^(1)12 



E 

K 



0(Pt - £x)^l{'K,ex ~ Mr) 



^kKJ 



dx 



J20{x-ex)rl{K,ex-x) 



K 



(66) 

X — fl-^ 



dx ^ ex 

k 



d{x - efc)Fo(k, ek - x) i 

ek — i?(K — k) — ^t + + X \x=ii^ 



However this is not yet the result we are looking for. Indeed this variation is at fixed /i^ whereas, just as in the 
preceding parts, we want to work in the canonical ensemble and find dpu^jdn^ at fixed n^, while Eq.(66l gives us 
d^l/dn\^ at fixed /i^. These two quantities are related by: 



( ^ 



From Eq.(56l and Eq.(60l we have: 



f 9nt \ 
\dnj 



kK 



9{n-f - ex)0{ek - Mt) 



Vl{K,ex-^i^) 



-iV;^Fo(A:;^,0) 



(67) 



(68) 



while, since we are in the limit of a vanishing number of polarons 
order result Eq.(61) for the polaron chemical potential: 



ZNf 



duP, , , 

Z-^ =iV;^Fo(fcF,0) + V 

d^i^ ^ 



0, d^i/dn^ is merely obtained from the lowest 
^(Mt - ^K)5{ek - /i^)F^(K, ex - Mr) 



kK 



E, 



(1) 

kK 



(69) 
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where Np = dfi^/dn^ = J^k ^(^fe ~ A*t) density of states at the spin-up Fermi level, and the factor 1/Z takes 

again into account that Fq depends on /z^. Our result is obtained from Eq.(66),Eq.(67),Eq.(68) and Eq.(69). 

It is interesting to consider the weak coupling limit which provides a consistency check oTour complicated result. 
Indeed in this case we have Fq = 2'Ka/mr and, as we have seen in Eq.(57l, S]-|-(k,w) is a constant which gives just a 
shift of the spin-up chemical potential. We have merely to replace fi^ by fi^ — 27r an^/m^. Otherwise the situation is 
not modified for the spin-up particles which behave still as non-interacting particles. Hence the situation is unchanged 
with respect to the calculation Eq.(61) to lowest order in n^. However, in the canonical ensemble, we are interested 
in the dependence of the polaron chemical potential fj.^ on the spin- up particle number n-^. Accordingly the formula 
giving in terms of n^, i.e. in terms of l/kpa, is unchanged in this weak coupling regime, and we should not find 
any correction. 

Indeed in this case we have [3] fi^ 



which is of order 



and dfi^/dn-^ = Fq in agreement with Eq.(69), as it should be, 
since from Eq.(65), we have Z ~ 1, the correction being of second order in Fq, i.e. in a^. In Eq.(68), the first term, 
is neghgible compared to the second one of order a, leading to dn^/dn^ 

' ' ^2 



a 



-NpTo. Hence the 

second term in the right-hand side of Eq.(67) is A^_fFq. On the other hand the first and third terms in the right-hand 
side of Eq.(66l are proportional to Fg, i.e. to a^. Hence in the weak coupling limit they are negligible compared to 



the second one which is proportional to Fq. Again we have only to take in this term the derivative of the Heaviside 
function and we find in this limit dfi^/dn^ 



Mi 



-NpTq. As a result, when we carry this value into Eq.(67), we find 



that dfii/dni 



, which is formally of order in our calculation, is actually exactly zero, as it should be. 



Similar simplifications are also present in the general case, and we can write the final result explicitly as: 
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(70) 



Here we have used our implicit notations k > kp and q < kp, q' < kp, while p runs over all possible values of the 
wavevector. We have set ojq — eq — /i-f , and Wpq = — ^^(p — q) + Finally, in the last term, the brackets denote 
the angular average on the direction of k running on the Fermi surface — fi^, which is easily performed analytically. 

Note that the apparent singularities in the first two terms, which would occur for -Eq/q = 0, are actually not present 



because the numerators of these first two terms cancel at the same place. 

Once the right-hand side of Eq.(70l is evaluated, the corresponding contribution E^^l to the interaction energy is 
given by: 



E, 



(3) 



^ Im. ( ^ 

n^Ep 2 Ep \dni 



(71) 



The coefficient in front of in Eq.(71) is plotted in Fig ji] for equal masses — to^. It is worth noting that the 
result is positive, a result expected physically but by no means obvious from the explicit expression Eq.(70). It should 
be noted that the result is quantitatively fairly small. For example at unitarity the coefficient is 0.041. This can be 
understood by noticing that the result goes very rapidly to zero in the weak coupling limit since we have seen that it 
behaves as in this limit. 



IV. DISCUSSION AND CONCLUSION 



In the preceding sections we have calculated contributions to the expansion in powers of the density x of the total 
system energy. They correspond to the x^ term and the x^/'^ term in this expansion. We note that this expansion 
is not restricted to integral powers x" as one would expect classically, but rather to integral powers x'^/^ of x^^^, 
that is we have actually an expansion in integral powers of the polaron Fermi sea radius pp. Naturally a few terms 
are missing at the beginning of this expansion and the first terms are the x term, coming from the polaron binding 
energy, and the x^^^ term, coming from the polaron kinetic energy. We have found terms proportional to x^ and x'''^^, 
and it is quite likely that all the higher powers of x-^^^ have nonzero coefficients. For example the polaron effective 
mass depends also on the polaron density, and we expect the first correction to the a; = result to be proportional to 
X (presumably corresponding to a decrease of the effective mass). This will lead to a x^^^ term (presumably with a 
positive coefficient). 
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FIG. 4: (Color online) E^'^ / (n^Epx^) as a function of l/kpa in the case of equal masses m-|- = mj,. 



It should be noted that we have not proceeded to a systematic expansion in powers of x within a single coherent 
framework. Proceeding in such a way would require a diagrammatic analysis. This is certainly a desirable goal, 
but unfortunately this looks fairly complicated to perform effectively. Our position has rather been to start from 
a physical point of view, and then to calculate each contribution in the more convenient way, either by the simple 
hamiltonian approach or diagrammatically. Nevertheless it is clear that, in such a diagrammatic analysis, the term 
calculated in section III B will appear as a polaron-polaron interaction term mediated by the spin- up Fermi sea (while 
the one calculated in section III A is just a direct manifestation of the polaron Fermi sea). Our lack of systematic 



expansion leaves open in principle the possibility that, up to the order we have considered, some terms are missing, 
although we believe that this is quite unlikely. 

We can now gather our results to get the total energy £ out of our results. We have to add to our interaction energy 
terms the energy of the free spin-up Fermi sea, the contribution from the chemical potential of isolated polarons and 
the kinetic energy of the non interacting polaron Fermi sea (2j . Focusing on the specific case of unitarity as in Ref. [5] , 
we have: 



0.6156 2:- 



3 X^/^ E; 



5 1.20 



(1) 



+ 



E: 



(2) 



E: 



(3) 



n^Ep 



(72) 



where the polaron binding energy and effective mass are taken from [7^ . The result is displayed in Fig. [s] and it is in 
clear disagreement with the Monte-Carlo results and the analytical curve of Ref. [2] as soon as x ^ 0.1. Naturally in 
this range it is reasonable to believe the numbers given by Monte-Carlo calculations since they take full account of 
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FIG. 5: (Color online) Total reduced energy £/{n^Ep) at unitarity as a function of the relative spin-down population x = ni/n^ 
(note that in Ref.[2] this is 5£ / {SnfEp) which is plotted). The (blue) thick dashed l ine is the expansion in powers of x limited 
to the x^^^ term with the coefficients from Ref.[7], i.e. th e th ree first terms of Eq.(72l. The (red) full line is this expansion 
going now up to the x"^^^ term, i.e. our full result Eq.(72|. The (black) dotted-dashed line is the same expansion, but 
omitting the dominant contribution coming from E^^^, and retaining only E^^^ and E^'^^. The (green) thin dashed line is the 
phenomenological formula Eq.(76l with p — 1.2. The blue dots are the Monte-Carlo results from Ref. [2]. 



all the effects coining from the polaron finite density (although they work with a small number of particles while in 
contrast our calculations are in the thermodynamical limit). 

The most natural way to understand this discrepancy is to remark that there is no reason to believe that, in our 
powers expansion, terms of higher order than the ones we have considered (i.e. going beyond x'^^^) do not contribute 

in a very important way. This is supported by the fact that the contribution from E^^^ is quite important. This finding 
implies that the series is not rapidly converging and that accordingly terms beyond this one give also an important 
contribution. Actually, at the start of the calculation, the only hope of agreement between our expansion and Monte- 
Carlo results was to find the opposite situation, that is all the terms we would find would be small. This would imply 
a rapidly converging series expansion which could be stopped at the order we have reached. Nevertheless this would 
have meant that polarons are weakly interacting objects, which is hard to understand physically since obviously the 
polarization cloud is easily perturbed. Hence, although our result is disappointing for quantitative agreement, it is 
quite reasonable physically. In this respect we note that, in Ref. [2], the reasonable agreement between Monte-Carlo 
results and the analytical result corresponding to the first three terms in Eq.( 72 ) was to a large extent due to the use of 
the value m* /m^ ~ 1 (obtained from Monte-Carlo calculations) for the polaron effective mass. Surprisingly the result 
is quite sensitive to the value of to*, and making use of m* /m^ = 1.20 leads to an analytical result fairly different from 
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the Monte-Carlo results. This sizeable difference can quite naturally be interpreted as the important overall effect of 



the interaction between polarons. We note incidentally that, if we omit the dominant interaction contribution i?j 



(1) 



and retain only El^^ and E^^l, we obtain a result in very good agreement with Monte-Carlo calculations as it can be 
seen on Figj5] 

There is naturally another possible source for the discrepancy we have found. Our calculation of the various 
coefficients is only approximate since we have stayed within the single particle-hole approximation. For the calculation 
of the polaron binding energy and effective mass, this approximation has proved to be remarkably accurate because 
taking more particle-hole excitations leads to a series which converges extremely rapidly toward the exact result [7]. 
It is likely that a similar conclusion applies for the coefficients we have calculated. This was our reason, in addition to 
simplicity, to stay at this level of approximation. However there is no good reason to believe that, if convergence there 
is, it is as fast as for the binding energy. The effective mass displays already a somewhat slower convergence. Hence it 
is quite possible that the convergence is not as fast for our coefffcients. In particular we have noticed, when calculating 
^int ' that keeping only a single particle-hole excitation was leading to quite inconvenient situations. Even if we have 
found our way around these problems, it seems likely they indicate that the single particle-hole approximation is not 
so accurate. Hence it is quite possible that the exact results for the power expansion coefficients are somewhat away 
from our findings. Nevertheless we believe it is quite unlikely that the corresponding changes reduce the interaction 
effects to such an extent that the polaron-polaron interaction could be considered as small enough to be neglected. 

We consider now the work of Mora and Chevy in the light of our results. First of all they have neglected the 
direct interaction term we have calculated in section |ll] since it is a x*"/^ term, which is a coherent point of view. 
Nevertheless, as we have stressed, this term becomes negligible only for very small x, in a range which is below the 
one involved in experiments or in Monte-Carlo calculations. We are then left with the term. Here we identify the 
term they have found with the term, imposed by thermodynamics, to go from the grand canonical to the canonical 
ensemble, i.e. the second term in the right-hand side of Eq.(67). Indeed, just as in Eq.(71), this term gives to the 
total energy a contribution: 



,(3) 



1 



^thc 



(73) 



Since, as we have seen, thermodynamics implies 9n^/97i^|^ = — 1^ , and diJi/drif 



(2/i-|-/3n-|-)9/i^/9/i^|^^, we obtain: 



^thcrm 



(74) 



which is just their result for the interaction term. This is coherent with their assumption that in the grand canonical 
ensemble, one has a mixture of two ideal Fermi gases of polarons and majority atoms, which means that there is 
no term describing their interaction. Hence in the canonical ensemble only the term required by thermodynamics 



appears. By contrast we have found a contribution in the grand canonical ensemble, namely the first term in Eq.(67l, 



and in the weak coupling limit a — >■ we find an overall result which behaves like in contrast with their a" 
behaviour. Finally they ascribe their result to the effect of Pauli blocking by the polaron Fermi sea wh ich they find 
behaving as x^. It is tempting to identify this effect with the one we have calculated in section III A but we have 



found a a;^/^ behaviour which is negligible in their framework. Finally we have mentioned that our"? term can be 
seen as an indirect interaction between spin-down atoms mediated by the spin-up Fermi sea. This is in agreement 
with the interpretation proposed by Yu, Zollner and Pethick [TB]. However our microscopic result Eq.(70) is much 



more complex than the one they propose within Fermi liquid theory. This might be due to the fact that, in this way, 
they restrict themselves to indirect interaction through low frequency perturbations, while we have not assumed such 
a restriction. It is also puzzling that they end up with the result obtained by Mora and Chevy, while we have seen 
that it has a purely thcrmodynamical interpretation. 

Finally it is interesting to present a simple phenomenological model which displays explicitly the troubles one meets 
when trying to perform a power expansion. In this model we treat at first symmetrically the spin-up and spin-down 
atoms. First there is the standard kinetic energy associated with the respective Fermi seas. Concentrating on the 
unitary case (and assuming also equal masses), we can in a simple phenomenology disregard the difference between the 
polaron effective mass and the bare mass. Introducing the total density n = n-^ + ni and the associated Fermi energy 
Ep = (67r^n)^/'^/2TO, the sum of the kinetic energies of the two Fermi seas is nEp (3/5)(m^/^ -I- w^/^) where u — n^/n 
and V = n^jn (and = 1). With respect to the polaron binding energy we know its value for — > (or — )• 0). 
However the effective binding energy decreases when increases because the polarization of the spin-up Fermi sea 
has to be shared between all the spin-down atoms (this is another way to see the polaron-polaron interaction). On 
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the other hand this binding energy will clearly be zero when = n since there is no spin-up Fermi sea anymore. 
Phenomenologically we may choose an interpolating function between these two limits, and write for example the 
interaction energy as —n^Ep eh(u) = —nEp u^/^vei,{v) with eb{v) — 0.6 (1 ~ vY, allowing for one fitting parameter p. 
Writing the corresponding expression for the spin-up polarons, we end up with the following expression for the total 
energy £ 

nEp 5 

The corresponding expression for Ejn^Ep in terms of a; = ^ih^ obtained by dividing the above expression by 
M^/^, and using u = 1/(1 + x) and v = x/{l + x). This leads to: 

It turns out that p = 1.2 gives a very good fit to the Monte-Carlo results as it can be seen from Fig[5] On the other 
hand since the radius of convergence for the series expansion of 1/(1 -|- a;) in powers of a: is 1, the expansion of £ in 
powers of x diverges for x = 1. Hence there is no way with a power expansion to have an agreement with Monte-Carlo 
results for a; = 1. And accordingly this cast some doubts on the possibility to get good results for lower values of x 
by making a series expansion and retaining a fairly large number of terms. 
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